
forvalues i=1910/1916 {
use "data\influenza_town_panel.dta",replace
local j=`i'+1
gen postx=(year>=`i'& year<=`j')

quietly tab did, gen(d_did)
forvalues x=1/205 {
gen d_did`x'_post=postx*d_did`x'
}

quietly reghdfe log_deaths_town d_did1_post-d_did205_post, absorb(cid) noconst
predict p
collapse (mean) p if year==`i' , by(indian_dm1918 did)
gen sample=1
save "temp_cdf_`i'.dta",replace
}

use "data\influenza_town_panel.dta",replace
tab did, gen(d_did)
forvalues i=1/205 {
gen d_did`i'_post=post*d_did`i'
}

quietly reghdfe log_deaths_town d_did1_post-d_did205_post, absorb( cid) noconst
predict p
collapse (mean) p if year==1918 , by(indian_dm1918 did sid )

gen s=1
forvalues i=1910/1916 {
append using "temp_cdf_`i'.dta"
}
replace indian_dm1918=4 if s==. & indian_dm1918==1
replace indian_dm1918=3 if s==. & indian_dm1918==0

drop if indian_dm1918==.

cdfplot p, by(indian_dm1918) graphregion(color(white)) xtitle("Bi-annual deviations from district average deaths") legend(label(1 "British DO 1918") label(2 "Indian DO 1918") label(3 "British DO 1910-17") label(4 "Indian DO 1910-17")) 
graph export "output\figure1.pdf", as(pdf) replace

forvalues i=1910/1916 {
erase "temp_cdf_`i'.dta"
}